Assessing injury on GOM cetaceans

1 Preamble

This is an Electronic Supplement to the manuscript Marques et al. “Quantifying Deepwater Horizon oil spill induced injury on pelagic cetaceans” submitted to Marine Ecology Progress Series (MEPS).

There are 8 Electronic Supplements to the paper. The master file containing links to all the other 7 additional Electronic Supplements related to this paper is ES0_ElectronicSupplements.

You might be reading this file as a pdf or as an html. The links on this file only work if you are using the html version of it, available via the github repository or if you compiled it yourself as html and you have all the 8 html files in the same folder. Otherwise, as a pdf distributed as an Electronic Supplement to the MEPS paper, the links might not work. They might work. If it is possible, we can work with the MEPS Editorial Office such that we can add links below that will link to actual files, say the pdfs of each of these 8 files, on the publisher server.

1.1 Version history

This section details the version history for static pdf files submitted as Electronic Supplement pdfs:

• 1.0 [12 Aug 2022] Version included as a pdf Electronic Supplement in the MEPS original submission

• 2.0 [10 Feb 2023] Version included as a pdf Electronic Supplement in the MEPS re-submission after 1st round of reviewer’s comments

2 Introduction

In this document we present the results of the age, sex and class structured model implemented for all species/stock in the GOM.

We report here the injury metrics originally considered in Schwacke et al. (2017): (1) lost cetacean years (LCY), the difference between the baseline and injured population sizes, summed over the entire modeled time period; (2) years to recovery (YTR), the number of years required before the injured population trajectory reaches 95% of the baseline population trajectory; and (3) maximum proportional decrease (MPD), the difference between the two population trajectories when the injured trajectory is at its lowest point, divided by the baseline. Note that LCY is intuitively the metric that is most dependent on initial population size.

3 Running the simulations

The code below runs the simulations for all species, and stores the results in appropriate files, stored in appropriate species specific folders.

This code is included for completeness and does not need to be run as part of compiling this dynamic report as the results are already available in relevant files and it takes a considerably log time to run.

Creating objects to hold summary results:

The taxonomic units and the corresponding codes considered in this document are:

• Bwsp - beaked whales, Beaked whales spp
• Fatt - pygmy killer whale, Feresa attenuata
• Ggri - Risso’s dolphin, Grampus griseus
• Gmac - short-finned pilot whale, Globicephala macrorhynchus
• Kosp - Kogia species, Kogia sp.
• Pele - melon-headed whale, Peponocephala electra
• Pmac - sperm whale, Physeter macrocephalus
• Satt - pantropical spotted dolphin, Stenella attenuata
• Sbre - rough-toothed dolphin, Steno bredanensis
• Scly - Clymene dolphin, Stenella clymene
• Scoe - striped dolphin, Stenella coeruleoalba
• Sfro - Atlantic spotted dolphin, Stenella frontalis
• Slon - spinner dolphin, Stenella longirostris
• Ttro - offshore bottlenose dolphins, Tursiops truncatus
• Ttrs - shelf bottlenose dolphins, Tursiops truncatus

We re-order the results tables upfront so that results are organized by alphabetical order of the 4 letter code used to describe the taxonomic unit (this is also the order that the species are reported in tables 1 to 3 of the offshore paper).

The data for initial population sizes and proportion exposed are imported from a txt file produced inside ES2_ElectronicSupplements.

4 Results by species

We represent in turn, for each species, the simulated population trajectories under the oil spill and under the baseline scenarios, as well as the histograms of the distributions of the 3 injury metrics for each species. To each histogram four dashed vertical lines were added:

• green - median value for the injury metric
• black - mean value for the injury metric
• blue - 95% confidence interval for the injury metric

4.1 Physeter macrocephalus (Pmac, sperm whale)

Population trajectories

Mean for each of the injury metrics

##         LCY         YTR         MPD 
## 2549.139215    9.725500    6.375102

4.2 Kogia sp. (Kosp, kogia)

Population trajectories

Mean for each of the injury metrics

##         LCY         YTR         MPD 
## 1429.418506    4.898900    5.185119

4.3 Stenella clymene (Scly, Clymene dolphin)

Population trajectories

Mean for each of the injury metrics

##         LCY         YTR         MPD 
## 1901.525708    0.000000    1.801292

4.4 Peponocephala electra (Pele, melon-headed whale)

Population trajectories

Mean for each of the injury metrics

##         LCY         YTR         MPD 
## 3071.523882    2.137900    4.168256

4.5 Tursiops truncatus (Ttro, offshore bottlenose dolphin)

Population trajectories

Mean for each of the injury metrics

##          LCY          YTR          MPD 
## 11369.564397     6.586900     5.658282

4.6 Globicephala macrorhynchus (Gmac, short finned pilot whale)

Population trajectories

Mean for each of the injury metrics

##         LCY         YTR         MPD 
## 1057.863912    1.156400    3.538044

4.7 Feresa attenuata (Fatt, pygmy killer whale)

Population trajectories

Mean for each of the injury metrics

##         LCY         YTR         MPD 
## 1249.503621    2.308200    4.194471

4.8 Grampus griseus (Ggri, Risso’s dolphin)

Population trajectories

Mean for each of the injury metrics

##         LCY         YTR         MPD 
## 2449.136975    7.290900    5.735186

4.9 Steno bredanensis (Sbre, rough-toothed dolphin)

Population trajectories

Mean for each of the injury metrics

##         LCY         YTR         MPD 
## 2434.417513    1.640900    4.006926

4.10 Stenella attenuata (Satt, pantropical spotted dolphin)

Population trajectories

Mean for each of the injury metrics

##          LCY          YTR          MPD 
## 33809.632317     0.635700     3.474522

4.11 Tursiops truncatus (Ttrs, shelf bottlenose dolphin)

Population trajectories

Mean for each of the injury metrics

##          LCY          YTR          MPD 
## 35190.477524     2.435900     4.191836

4.12 Stenella longirostris (Slon, spinner dolphin)

Population trajectories

Mean for each of the injury metrics

##          LCY          YTR          MPD 
## 17714.767745    11.196700     8.899829

4.13 Stenella coeruleoalba (Scoe, striped dolphin)

Population trajectories

Mean for each of the injury metrics

##         LCY         YTR         MPD 
## 3715.207476    7.108100    5.883191

4.14 Stenella frontalis (Sfro, Atlantic spotted dolphin)

Population trajectories

Mean for each of the injury metrics

##         LCY         YTR         MPD 
## 7509.270139    0.000000    1.296869

4.15 Beaked whales spp (Bwsp)

Population trajectories

Mean for each of the injury metrics

##         LCY         YTR         MPD 
## 1320.899319    1.592200    3.816632

5 Comparisons across species

In this section we present some summary analysis of the results across the 15 taxonomic units considered.

5.1 Initial population and proportion exposed

The population size and the proportion exposed are \(a priori\) expected to be key determinants of injury, since the population size will have a direct impact on LCY and the proportion exposed should have a direct impact on all 3 measures.

Initial population size (N) and proportion exposed (Pexp) for the different stocks considered.

Sp	Nstart	Pexp
Bwsp	3098	0.140
Fatt	2152	0.145
Ggri	3063	0.197
Gmac	2065	0.120
Kosp	2322	0.197
Pele	5784	0.152
Pmac	2561	0.191
Satt	81233	0.155
Sbre	4867	0.148
Scly	9065	0.082
Scoe	5011	0.247
Sfro	48688	0.058
Slon	16501	0.377
Ttro	15791	0.206
Ttrs	64897	0.177

The proportion exposed against the initial population size is shown here for the 15 taxonomic units

5.2 Injury metrics

In the following table we present the means of the 3 injury metrics:

Injury results for the different stocks considered.

Sp	LCY	MPD	YTR
Bwsp	1321	-0.038	1.6
Fatt	1250	-0.042	2.3
Ggri	2449	-0.057	7.3
Gmac	1058	-0.035	1.2
Kosp	1429	-0.052	4.9
Pele	3072	-0.042	2.1
Pmac	2549	-0.064	9.7
Satt	33810	-0.035	0.6
Sbre	2434	-0.040	1.6
Scly	1902	-0.018	0.0
Scoe	3715	-0.059	7.1
Sfro	7509	-0.013	0.0
Slon	17715	-0.089	11.2
Ttro	11370	-0.057	6.6
Ttrs	35190	-0.042	2.4

In the following table we present the medians and the 95% confidence intervals for our current estimates of injury. These correspond to the results in table 3 in the main paper.

Median injury results for the different stocks considered and respective 95% confidence intervals by the percentile method

Sp	LCYmed	LCYl	LCYu	MPDmed	MPDl	MPDu	YTRmed	YTRl	YTRu
Bwsp	1207	445	2899	3.683	1.386	6.951	0	0	11
Fatt	1122	437	2809	4.186	1.939	6.642	0	0	13
Ggri	2264	978	5049	5.623	2.569	9.494	8	0	18
Gmac	964	385	2291	3.435	1.486	6.216	0	0	12
Kosp	1294	519	3112	5.114	2.315	8.492	5	0	14
Pele	2761	1081	6767	4.140	1.925	6.601	0	0	12
Pmac	2356	996	5240	6.344	2.672	10.189	11	0	21
Satt	31372	12884	67606	3.433	1.563	5.615	0	0	9
Sbre	2232	901	5207	3.995	1.836	6.349	0	0	11
Scly	1726	664	4067	1.734	0.747	3.261	0	0	0
Scoe	3387	1369	7897	5.673	2.479	10.426	8	0	17
Sfro	6961	2722	15509	1.264	0.542	2.280	0	0	0
Slon	15255	5302	44903	8.363	3.277	17.716	12	0	21
Ttro	10537	4597	23220	5.610	2.546	9.125	8	0	16
Ttrs	32584	13377	71967	4.152	1.859	6.758	0	0	13

These are the values that are used in table 3 of the paper, and the code to do so is in the .Rmd.

Representing the above information visually, focusing on each of the 3 injury metrics:

6 Exploring the injury metric results

In this section we explore the results obtained in terms of the different injury metrics and how these injury metrics might be explained by the different input parameters. To do so we consider linear models to model the mean value of the injury per species as a function of the mean values of the input parameters. After a preliminary exploratory analysis, only the following four parameters seem to have a non-negligible potential relevant effect on the injury metrics:

• Proportion Exposed
• Initial population size
• Gestation duration
• Survival reduction

6.1 Exploratory analysis

Here we plot the values of the injury metrics as a function of the four key parameters going into the model.

6.1.1 Proportion exposed

The next three plots describe the relationship between the proportion exposed and Lost Cetacean Years, the Maximum Proportion Decrease and the Years to recovery, respectively

6.1.2 Initial population size

The next three plots describe the relationship between the initial population size and Lost Cetacean Years, the Maximum Proportion Decrease and the Years to recovery, respectively.

6.1.3 Gestation duration

The next three plots describe the relationship between the gestation duration and Lost Cetacean Years, the Maximum Proportion Decrease and the Years to recovery, respectively

6.1.4 Survival reduction

The next three plots describe the relationship between the survival reduction and Lost Cetacean Years, the Maximum Proportion Decrease and the Years to recovery, respectively

6.2 Explaining injury metrics across species

Here we consider standard linear models to explain the injury metrics as a function of the input parameters. To obtain what we considered the most parsimonious model we conducted a traditional model selection approach using AIC, selecting the model with the lowest AIC.

We considered only the 4 explanatory variables which, given the exploratory data analysis above, could explain the injury metrics. Given only 4 explanatory variables, we conducted a full search over the possible models, ignoring interactions, using function bestglm from the package bestglm (McLeod & Lai, 2020).

6.2.1 LCY

The summary of the best linear model follows.

## AIC
## BICq equivalent for q in (0.00257603089392322, 0.715027970358707)
## Best Model:
##                  Estimate     Std. Error   t value         Pr(>|t|)
## (Intercept) -8585.8086023  2457.09226124 -3.494296 0.00442815340683
## N               0.4263961     0.03567448 11.952412 0.00000005054783
## p           54656.6056191 12278.22336920  4.451508 0.00079078410382

## 
## Call:
## lm(formula = y ~ p + N, data = LCY4bestglm)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -7835.4 -1189.7   855.8  1477.7  6430.2 
## 
## Coefficients:
##                Estimate  Std. Error t value     Pr(>|t|)    
## (Intercept) -8585.80860  2457.09226  -3.494     0.004428 ** 
## p           54656.60562 12278.22337   4.452     0.000791 ***
## N               0.42640     0.03567  11.952 0.0000000505 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3380 on 12 degrees of freedom
## Multiple R-squared:  0.9263, Adjusted R-squared:  0.9141 
## F-statistic: 75.45 on 2 and 12 DF,  p-value: 0.0000001598

We can look at the ANOVA table for the most parsimonious model based on lowest AIC criteria.

## Analysis of Variance Table
## 
## Response: y
##           Df     Sum Sq    Mean Sq  F value        Pr(>F)    
## p          1   91857453   91857453   8.0395       0.01502 *  
## N          1 1632291730 1632291730 142.8602 0.00000005055 ***
## Residuals 12  137109617   11425801                           
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Using the most parsimonious model, the proportion of the variability explained by each input retained in the regression model and the remaining unexplained variation can be obtained from the ANOVA table of the linear regression (as shown above). The proportion of variance explained in Lost Cetacean Years by the initial population size is 87.7 %, while the proportion exposed accounts for 4.94 %. Only 7.37 of the variability remains unexplained.

We illustrate how well the model predicts the observed values in the plot below:

6.2.2 YTR

The summary of the best linear model follows.

## AIC
## BICq equivalent for q in (0.0565670801885912, 0.714907105110615)
## Best Model:
##                 Estimate  Std. Error   t value       Pr(>|t|)
## (Intercept) -14.46579200 3.864204450 -3.743537 0.002803820209
## gd            0.02778531 0.009303748  2.986464 0.011348598337
## p            44.99620475 5.546928494  8.111914 0.000003260501

## 
## Call:
## lm(formula = y ~ p + gd, data = YTR4bestglm)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -2.364 -1.229  0.256  1.220  2.177 
## 
## Coefficients:
##               Estimate Std. Error t value   Pr(>|t|)    
## (Intercept) -14.465792   3.864204  -3.744     0.0028 ** 
## p            44.996205   5.546928   8.112 0.00000326 ***
## gd            0.027785   0.009304   2.986     0.0113 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.515 on 12 degrees of freedom
## Multiple R-squared:  0.8502, Adjusted R-squared:  0.8253 
## F-statistic: 34.06 on 2 and 12 DF,  p-value: 0.00001129

We can look at the ANOVA table for the most parsimonious model based on lowest AIC criteria.

## Analysis of Variance Table
## 
## Response: y
##           Df  Sum Sq Mean Sq F value      Pr(>F)    
## p          1 135.955 135.955  59.199 0.000005592 ***
## gd         1  20.483  20.483   8.919     0.01135 *  
## Residuals 12  27.559   2.297                        
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

The proportion of variance explained in Years to Recovery by the proportion exposed is 73.89 %, while the gestation duration accounts for 11.13 %. About 14.98 % of the total variability remains unexplained.

We illustrate how well the model predicts the observed values in the plot below:

6.2.3 MPD

The summary of the best linear model follows.

## AIC
## BICq equivalent for q in (0.000211081585492212, 0.79448839248248)
## Best Model:
##                   Estimate    Std. Error    t value              Pr(>|t|)
## (Intercept) -0.14127227832 0.03051310684  -4.629888 0.0007285396508132729
## gd          -0.00006689489 0.00001223598  -5.467064 0.0001957224385135165
## sr           0.19516442031 0.03283771043   5.943302 0.0000968549623446209
## p           -0.25059730715 0.00616856884 -40.624870 0.0000000000002441672

## 
## Call:
## lm(formula = y ~ p + gd + sr, data = MPD4bestglm)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -0.0020501 -0.0011974  0.0001799  0.0007734  0.0033525 
## 
## Coefficients:
##                Estimate  Std. Error t value          Pr(>|t|)    
## (Intercept) -0.14127228  0.03051311  -4.630          0.000729 ***
## p           -0.25059731  0.00616857 -40.625 0.000000000000244 ***
## gd          -0.00006689  0.00001224  -5.467          0.000196 ***
## sr           0.19516442  0.03283771   5.943 0.000096854962345 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.001685 on 11 degrees of freedom
## Multiple R-squared:  0.9936, Adjusted R-squared:  0.9918 
## F-statistic: 567.1 on 3 and 11 DF,  p-value: 0.000000000002469

We can look at the ANOVA table for the most parsimonious model based on lowest AIC criteria.

## Analysis of Variance Table
## 
## Response: y
##           Df    Sum Sq   Mean Sq  F value             Pr(>F)    
## p          1 0.0044337 0.0044337 1561.410 0.0000000000003305 ***
## gd         1 0.0002966 0.0002966  104.451 0.0000005944763062 ***
## sr         1 0.0001003 0.0001003   35.323 0.0000968549623446 ***
## Residuals 11 0.0000312 0.0000028                                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

The proportion of variance explained in MPD by the proportion exposed is 91.19, while the gestation duration and the survival reduction only account for only 6.1 % and 2.06 %, respectively. About 0.64 % of the variability remains unexplained.

We illustrate how well the model predicts the observed values in the plot below:

7 References

DWH MMIQT 2015, Models and analyses for the quantification of injury to Gulf of Mexico cetaceans from the Deepwater Horizon Oil Spill, MM_TR.01_Schwacke_Quantification.of.lnjury.to.GOM.Cetaceans LINK.

McLeod A, Xu C, Lai Y (2020). bestglm: Best Subset GLM and Regression Utilities. R package version 0.37.3, https://CRAN.R-project.org/package=bestglm.

Schwacke, L.H., L. Thomas, R.S. Wells, W.E. McFee, A.A. Hahn, K.D. Mullin, E.S. Zolman, B.M. Quigley, T.K. Rowles and J.H. Schwacke. 2017. An age-, sex- and class-structured population model for estimating nearshore common bottlenose dolphin injury following the Deepwater Horizon oil spill. Endangered Species Research 33: 265-279. DOI: 10.3354/esr00777.